SEd = input('What is the sample standard error? ');
meand = input('What is the sample mean? ');
n = input('How many subjects were run? ');
 
Vard = n*SEd^2;
SSd = Vard*(n-1);
 
likelihoodmax = 0;
theta(1) = meand - 5*SEd;
inc = SEd/100;
 
 for B = 1:1000
    theta(B) = theta(1) + (B-1)*inc;
    likelihood(B) = (SSd + n*(meand - theta(B))^2)^(-n/2);
   
    if likelihood(B) > likelihoodmax
        likelihoodmax = likelihood(B);
    end
 end
 for B = 1:1000
     likelihood(B) = likelihood(B)/likelihoodmax;
 end
 
 outofrange = meand - 6*SEd;
 
 begin8 = outofrange;
 begin32 = outofrange;
 end8 = outofrange;
 end32 = outofrange;
 for B = 1:1000
     if begin8 == outofrange
         if likelihood(B) > 1/8
             begin8 = theta(B);
         end
     end
     if begin32 == outofrange
         if likelihood(B) > 1/32
             begin32 = theta(B);
         end
     end
     if and(begin8 ~= outofrange, end8 == outofrange)
         if likelihood(B) < 1/8
             end8 = theta(B);
         end
     end
     if and(begin32 ~= outofrange, end32 == outofrange)
         if likelihood(B) < 1/32
             end32 = theta(B);
         end
     end
 end
 
 
  theta1 = input('What is the population mean assumed by the first hypothesis? ');
  theta2 = input('What is the population mean assumed by the second hypothesis? ');
  B1 = int16((theta1 - theta(1))/inc  + 1);
  B2 = int16((theta2 - theta(1))/inc  + 1);
  likelihoodratio  = likelihood(B1)/likelihood(B2)
    
     begin32
     end32
     begin8
     end8